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Abstract:  Katzenelson’s  algorithm  and  its  variants  are  powerful  tools 
for  solving  nonlinear  networks  which  are  modeled  by  piecewise-linear 
4PUD"- characteristics.  But,  when  nonlinear  network  sizes  become  very 
large  such  as  in  VLSI  chip  cases,  excessive  epu  time  and  storage  are 
required  during  the  solution  process  using  Katzenelson' s  algorithm. 
Decomposition  techniques  are  necessary  in  the  analysis  of  VLSI  cir¬ 
cuits.  Nonlinear  Gauss-Seidel  iterative  methods  are  often  adopted  in 
solving  large  decomposed  system  of  equations.  However,  Nonlinear 
Gauss-Seidel  iterative  process  will  converge  under  certain  condi¬ 
tions.  The  combination  of  Katzenelson  and  Gauss-Seidel  methods  pro¬ 
posed  here  takes  advantages  of  both  Katzenelson  and  the  Gauss-Seidel 
methods.  It  decomposes  the  whole  network  into  small  subcircuits  by 
Gauss-Seidel  method  and  solves  these  small  subcircuits  by 
Katzenelson's  algorithm  separately  (or  even  these  subcircuits  can  be 
solved  by  Katzenelson's  algorithm  at  same  time  with  parallel 
processors, if  Jacobi  Method  is  used  as  decomposition  technique).  The 
convergence  properties  of  the  method  is  studied  in  detail,  and  exam¬ 
ples  are  given  here  to  illustrate  the  approach  in  the  dc  analysis  of 
bipolar  and  MDS  transistors  circuits.  Jj.  ’)  '•  °  I  0/ 
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1  Introduction 

The  spectacular  growth  in  scale  of  integrated  circuit  in  the 
VLSI  era  needs  aore  efficient  methods  for  circuit  simulation.  In  the 
area  of  computer-aided  nonlinear  network  analysis,  the  techniques  of 
piecewise-linear  ( PWL)  approximation  and  analysis  are  efficient  in 
computation,  and  have  received  considerable  interest  during  the  last 
decade.  A  well-known  technique  namely,  Katzenelson' s  algorithm  [1] 
which  was  originally  applied  to  nonlinear  resistive  network  with 
two-terminal  monotonic  elements,  has  been  extended  to  solve  more  gen¬ 
eral  cases  [2], [3]. 

In  piecewise-linear  analysis,  the  nonlinear  element  charac¬ 
teristics  are  approximated  by  continuous  piecewice-1 inear  functions 
and  the  nonlinear  network  model  can  be  expressed  as 

F(X)  -  J(b)X  +  W<B>  -  Y  m**l , . . . , r  (1) 
where  F(X)  is  a  continuous  P1TL  mapping  from  Ra  into  itself,  X  is  a 
point  in  Rn  and  represents  a  set  of  chosen  variables  in  a  given  net¬ 
work  and  Y  is  an  arbitrary  point  in  Rn  and  represents  the  inputs,  the 
space  Rn  is  divided  into  a  finite  number  (r)  of  polyhedral  regions 
by  a  finite  number  of  hyperplanes,  is  a  constant  n  x  n  Jacobian 
matrix  and  is  a  constant  vector  defined  in  a  given  region  i,  and 
r  denotes  the  total  number  of  regions. 

The  PWL  equation  (1)  can  be  solved  by  one  of  several  methods: 
Newton's  method,  Katzenelson' s  algorithm,  Jacobi  iterative  method, 
Gauss-Seidel  and  its  extension  SOR  iterative  methods.  Newton's 
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method  suffers  nonconvergence  problems  unless  the  initial  guess  is 


close  enough  to  the  solution.  Katzenelson' s  algorithm,  which  is  s 
modified  algorithm  of  Newton's  method,  is  a  powerful  tool  for  solving 
continuous  PWL  equations  (1)  and  has.  good  convergence  properties  and 
has  been  refined  and  extended  to  more  general  cases.  It  has  been 
shown  that  as  long  as  all  the  Jacobian  matrix  determinants,  det 
..  ,r  in  (1),  have  the  case  sign,  there  exists  at  least  one 
solution  to  the  equation  (1)  and  Katzenelson's  algorithm  always  con¬ 
verges  to  a  solution  [31.  Even  the  sign  restriction  was  later  removed 
in  the  generalized  Katzenelson  method  [4]. 

The  basic  Katzenelson's  algorithm  [5]  is  as  follows: 

1.  Choose  an  initial  guess  X^,  k*0 

2.  Find  the  corresponding  region  k  and  construct  J^k*  and 

3.  Find  J(k)Xk  +  f(k>  «  Y(k) 

4.  Solve  J(k)AX(k)  -  T  -  Y(k)  «  AY(k) 

5.  Put  X(k+1)  »  X(k)  +  AX(k) 

6.  If  X<k+1)  is  in  region  k,  a  soluiotn  is  found;  otherwise, 
choose  X(k+1)  -  X(k>  +  \AX(k>  (0  <  \  £  i)  to  be 

on  the  boundary  of  region  k  and  cross  into  adjacent  region. 

7.  put  Y(k+1)  -  (1  -  \  )Y(k)  and  update  J(k)  to  J<k+1) 
using  a  dyad  relationship  and  return  to  step  4. 

However,  when  the  size  of  a  nonlinear  network  beoomes  very  large,  for 
instance,  when  the  number  of  nodes  n  >  ioo,  the  Jacobian  matrix  J^“* 
in  (1)  becomes  larger  than  100  x  100  and  the  number  of  regions  r  may 
be  also  very  large,  which  results  in  excessive  cpn  time  and  storage 
during  the  Katzenelson's  procedure. 
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Some  alternative  methods  for  solving  nonlinear  equations  are 
Gauss-Seidel.  Jacobi  and  SOR  iterations.  These  methods  partition  the 
n  dimensional  equation  into  n  one-dimensional  equations.  Each  equa¬ 
tion  is  then  solved  for  one  variable  by  assuming  the  other  (n-1) 
variables  known.  If  the  Gauss-Seidel  method  is  applied  to  solve  the 
equations,  equation  (1)  can  be  expressed  as  follows 

fl(xl““*  xi . xn}  "  *1 

•  •  • 

•  •  • 

£i(*1....»  *!.•••.  *a>  -  7i  (2) 


fn(x1,....  xj . xn)  -  yn 

where  f^,  and  yj  are  the  ith  component  of  function  F(X)  .vector  X 
and  vector  T  respectively.  The  Gauss-Seidel  iterative  formula  [6] 
becomes 


fi(x 


(k+1) 

1 


(k+1) 
xi-l  * 


xi 


(3) 


_(k+l) 

Xi 


“x, 


(4) 


When  solving  (3)  and  (4),  X^*(xi^  . . . . ,  xj*^ ,. . .  ,x^  are  known 
from  the  kth  iterate  and  xj^*^ ,. . . ,x are  also  known  from 
solving  the  first  (i-1)  components  during  the  (k+l)th  iterate.  Thus 
(3)  becomes  a  one-dimensional  equation  which  is  easy  to  solve  for 
Xj.  The  process  is  repeated  until  the  iterations  converge,  provided 
p(J^)  <  1,  where  p(J^)  is  spectral  radius  of  defined  as  the 
maximum  modulus  of  the  eigenvalues  of  matrix  [61. 


The  Gauss-Seidel  as  well  as  Jacobi  and  SOH  iterations  are 
guaranteed  to  converge  when  any  one  of  the  following  three  conditions 
is  satisfied  [8]:  (  These  conditins  are  sufficient  but  not  necessary) 

1.  J^  is  symmetric  and  positive  definite, 

2.  J(a,)  is  a  class  M-aatrix 

3.  J^  is  a  diagonally  dominant 

In  nonlinear  transistor  network  analysis,  the  second  and  third  condi¬ 
tions  are  more  applicable  sinoe  the  Jacobian  aatrix  is  generally  non- 
syaaetric.  A  class  M-aatrix  and  diagonally  dominant  are  defined  as 
follows  [6] ? 

Definition  1.1 

A  matrix  A  €  L(Rn)  is  an  M-aatrix  if  Awta^)  is  invertible, 

A  1*(b.j),  «ad  b-j  i.  0  for  all  i,  j*l,...,n  and  ajj  £  0  for  all 
i,j-l,...,n.  i*j. 

Definition  1.2 

A  matrix  A  €  L(Cn)  is  diagonally  doainant  if 

**44!  1  Z  l»iil  j-l,...,n 
j*i  1J 

and  strictly  diagonally  doainant  if  strict  inequality  holds. 

The  matrices  generated  by  the  nodal  approach  froa  networks  which 
consist  of  only  two- terminal  uncoupled  resistors  and  independent 
sources  are  class  M-aatrices  as  well  as  diagonally  dominant.  Vhen 
the  circuit  contains  also  linear  two-terminal  capacitors  and  if  the 
capacitors  are  replaced  by  companion  aodels  using  an  implicit 
integration  formula  its  nodal  aatrix  is  also  an  M-aatrix.  However, 
in  aany  transistor  networks,  definition  1.1  is  not  satisfied;  for 
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instance,  in  the  TIL  NAND  gate  shown  in  Fig.l.  In  order  to  fora  the 
Jacobian  matrix  for  the  TIL  NAND  gate  in  Fig.l,  we  use  the  Ebers-Moll 
model  shown  in  Fig. 2  to  replace  each  npn  bipolar  transistor  in  Fig.l. 
Each  diode  in  Fig. 2  is  modeled  by  piecewise-1 inear  characteristics 
shown  in  Fig. 3  (a)  so  that  the  PWL  model  for  npn  bipolar  transistor 
which  is  shown  in  Fig. 3  (b)  is  obtained.  The  sign  matrix  in  Fig. 4  is 
the  sign  of  elements  in  the  Jacobian  matrix  of  the  TIL  NAND  gate  in 
Fig.l,  where  +  means  the  element  in  this  position  conld  be  negative 
or  positive  depending  on  the  piecewise  region  in  which  the  nonlinear 
resistors  are  located;  snch  as,  the  node  voltages  Vj-4.8,  v^*4.9, 
▼5=5.0,  others~0.0,  the  signs  of  the  elements  in  row  4  and  colnmn  2, 
and  row  3  colnmn  4  are  positive,  so  the  Jacobian  matrix  is  not  a 
class  K-matrix.  The  Jacobian  matrices  are  not  always  diagonally  dom¬ 
inant  either,  since  controlled  sources  or  feedback  loops  may  exist  in 
some  circuits.  The  point  Gauss-Seidel  iteration  may  not  be  converge 
when  applied  to  solve  this  kind  of  circuits. 

Consider  a  nonlinear  network,  such  as  the  one  in  Fig. 5,  which 
consists  of  TIL  NAND  gates.  The  network  matrix  is  neither  an  M- 
matrix  nor  diagonally  dominant,  and  the  nonlinear  point  Gauss-Seidel 
iteration  does  not  converge  to  a  solution.  On  the  other  hand, 
Katzenelson' s  algorithm  when  applied  to  solve  this  network  converges 
to  the  solution.  But,  as  mentioned  earlier,  applying  Katzenelson' s 
algorithm  on  the  entire  circuit  would  require  a  relatively  excessive 
computational  time  and  storage.  In  some  cases  even  when  nonlinear 
point  Gauss-Seidel  method  converges  to  a  solution,  excessive  compute- 
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tional  tine  is  required  before  tbe  solution  is  reached.  In  this 
case,  a  block  Gauss-Seidel  iterative  method  could  converge  to  the 
solution  in  less  time. 

In  this  paper  a  PWL  block  Gauss-Seidel  approach  is  proposed.  The 
approach  combines  Katzenelson' s  algorithm  and  Gauss-Seidel  iteration 
for  solving  large  scale  nonlinear  networks  which  are  modeled  by  PfL 
functions.  For  example,  the  entire  network  in  Fig.S  is  first  decom¬ 
posed  into  small  subnetworks  (or,  subcircnits)  as  shown  in  Fig.l. 
The  subcircuits  are  solved  one  at  a  time  in  a  given  determined 
sequence.  The  aMthod  will  be  refered  to  as  Katzenelson-Gauss-Seidel 
method,  which  will  be  described  in  detail  in  section  3.  In  section 
2,  we  first  introduce  some  basic  properties  of  nonlinear  point 
Gauss-Seidel  iteration. 

2  Convergence  conditions  for  nonlinear  Point  Gauss-Seidel  iteration 

In  this  section,  we  introduce  some  definitions  and  theorems  and  give 
the  convergence  conditions  for  nonlinear  point  Gauss-Seidel  iteration 
followed  by  the  convergence  conditions  for  Katzenelson-Gauss-Seidel 
method.  When  considering  convergence  conditions,  certain  measures 
for  comparing  vectors  and  matrices  are  needed.  Comparing  vectors  in 
Rn  element  by  element  has  been  found  to  be  one  of  the  most  useful 
measures.  This  can  be  done  by  means  of  the  natural  (or  component¬ 
wise)  partial  ordering  on  Rn  defined  by 

For  X, Y  €  Rn,  X  <  I,  if  tnd  only  if  x,  i  y,  i-1 . n 


where  X,  T  ere  vectors  on  Rn  end  Xj,  y^  are  the  ith  components  of  X. 
Y  respectively.  Similarly,  for  nttrices 

For  A.  B  €  L(Rn,  R“) ,  A  <  B,  if  end  only  if  Sjj  i  hjj 

i”l». . . #■#  j“l,...n 

where  s^  end  h-j  ere  the  elements  in  row  i  and  colnan  j  in  matrices 
A  and  B.  From  above.  A  >  o  means  that  a^j  20,  and  A  is  called  a 
nonnegative  matrix.  We  use  F:  D  Rn->R®,  f^X),  i«l,...,m,  X  6  Rn»  as 

a  notation  for  a  mapping  F  with  domain  D  in  Rn  and  range  F(D)  in  R® 
and  with  components  f-,..4>f  .  The  following  are  some  definitions 
[6]  which  are  necessary  for  understanding  the  convergence  conditons 
of  the  Gauss-Seidel  iteration. 

Definition  2.1 

the  mapping  F:  DCRn  ->  R®  is  isotone  (  or  antitone)  (on  Rn)  if 

X  <  Y  for  any  X,  Y  €  D,  implies  that  F(X)  i  F(Y)  (or  F(X)  2  F(Y)  ); 
if  the  inequality  holds,  F  is  strictly  isotone  (or  antitone). 
Definition  2.2 

the  mapping  F:  DC  Rn  ->  R®  is  inverse  isotone  (on  Rn)  if 
F(X)  <  F(Y)  for  any  X,  Y  €  D,  implies  that  X  £  Y 

Definition  2.3 

the  mapping  F:  DCRn  ->  R®  is  diagonally  isotone  (on  Rn)  if  for 
any  X  €  Rn  the  functions 

*i,i:  CteR1IX+teie  D}->R1,  ^i>i(t)-Fi(X+tei).  i-l,...,n 

are  isotone,  where  e*  is  a  unit  basic  vector  with  the  ith  componont 
one  and  all  others  zero.  Analogously,  we  can  get 


Definition  2.4 


the  sapping  F:  DCS0  ->  Rn  ia  off-diagonally  antitone  (on  S°)  if 
for  any  I  6  Ea  the  functions 

^  (t€R1IX+te^€  D  }->R1,  ^  jUJ-FjU+te^)  i=j .  i.  j-1. . . .  .n 

are  antitone,  where  e^  is  a  unit  basic  vector  with  the  jth  componont 
oni  and  all  others  zero.  The  following  gives  the  definition  of  the 
M-function,  which  is  an  extension  of  the  M-matrix  to  the  nonlinear 
case. 

Definition  2.5 

the  mapping  F:  DCRn->Rn  is  an  M-function  if  F  is  inverse  iso¬ 
tone  and  off-diagonally  antitone. 

The  M-matrix  defined  in  section  1  has  the  properties  that  the  diago¬ 
nal  elements  are  positive  and  the  off-diagonal  elements  are  negative 
and  the  inverse  of  M-matrix  is  nonnegative.  The  M-function  has  simi¬ 
lar  properties.  The  following  theorem  describes  some  properties  of 
M-functions. 

Theorem  2.1 

Let  F:  DCR°  ->  Rn  be  an  M-function  (and  hence  injection)  .then 
F  and  F-*:  F(D)CR°  ->  Rn  are  strictly  diagonally  isotone.  If  F:  R°-> 
R°  is  continuous  and  surjective,  then  F  and  F“*:  Rn  ->  R°  are  surjec¬ 
tively  diagonally  isotone. 

The  proof  can  be  found  in  [6]. 
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With  these  definitions  and  theorem,  we  can  state  the  convergence 
properties  of  nonlinear  Ganss-Seidel  iteration. 

Theorem  2.2 

Let  F:  DC  Rn->Ra  be  continnons  of f-diagonally  antitone  and 
strictly  diagonally  isotone.  Suppose  that  for  some  B  €  Rn  there  exist 
points  X°,  Y°  €  Rn  such  that 

X°  <  T°,  FU°)  <  B  £  F(Y°> 

Then,  the  Gauss-Seidel  iterates  {Yk}  {X^}  given  by  (3)  and  (4)  and 
starting  from  Y°  and  X°,  are  uniquely  defined  and  satisfy 

X°  <  Xk  <  Xk+1  <  Yk+1  <  Yk  <  y°,  F(X)k  <  B  *  F(Y)k 

k*0 , 1 , . . .  (5) 

as  well  as 

limXk  =  X*  <  Y*  =limYk,  F(X*)  =  F(Y*)  =  B. 
k->«  k->» 

This  proof  can  be  found  in  [6],  and  is  given  in  the  Appendix  for  easy 
reference. 

Theorem  2 .3 

Let  F:  Rn->Rn  be  a  continuous  M-function  from  Rn  onto  itself. 
Then  for  any  B  6  Rn,  any  starting  point  X°,  the  Gauss~Seidel  iterates 
(3)  converge  to  the  unique  solution  X*  of  F(X)=B. 

The  proof  is  also  given  in  the  Appendix.  What  has  been  discussed  so 
far  is  suitable  for  general  continuous  nonlinear  functions,  which 
include  continuous  PWL  functions.  The  following  theorem  deals  with 


continuous  PWL  mappings: 


Theorem  2.4 


Let  F:  Rn->Rn  be  *  continuous  piecewice-linear  affine  mapping 
(1),  for  which  all  matrices  {J}  are  M-matrices,  then  F  it  an  M- 
f unction. 

This  theorem  is  proved  in  [10]. 

By  theorems  2.3  and  2.4.  as  long  as  the  matrices  (J)  in  piece- 
wise  mapping  (1)  are  M-matrices,  nonlinear  Gauss-Seidel  iteration  can 
be  used  in  PWL  cases  and  converges  to  the  unique  solution. 

Another  condition  which  guarantees  the  convergence  of  point 
Gauss-Seidel  iteration  is  diagonal  dominance.  If  a  matrix  A  is  diag¬ 
onally  dominant,  and  A  can  be  split  as: 

A  -  D  -  L  -  U 

where  D,  L  and  U  are  the  diagonal  part,  lower  triangular  part  and 
upper  triangular  part  of  A  respectively,  then  the  sequence  generated 
by  the  point  Gauss-seidel  iterative  method  in  solving  the  linear 
matrix  equation 

AX  -  P 

can  be  written  as: 

Xk+1  -  (D  -  L)_1UXk  +  (D  -  L)_1P  (6) 

provided  (D  -  L)  is  nonsingular. 

The  point  Gauss-Seidel  iterative  sequence  {Xk}  in  (6)  will  converge 
to  the  unique  solution  by  the  Diagonal  Dominance  Theorem  in  [7]  since 
the  matrix  A  is  diagonally  dominant. 

For  piecewise  linear  cases,  if  the  matrices  in  (1)  are 


diagonally  dominant  and  can  be  split  into: 
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j(m)  m  p(m)  _  ^(b)  _  n(m) 

Let 

M(m)  -  (D(b)  - 

and 

p(«),(D(«>  _  i/m) )_1(T  -  f<*>) 

provided  (D^  -  L^*^)  it  nonaingnlar.  If  the  Ganas-Seidel  iterative 

method  is  used,  equation  (1)  can  be  rewritten  as  follows: 

xk  „  M(m)xk-1  +  p(m)  (7) 

where  and  P^B^  have  different  valnes  depending  on  the  variable 

values  in  X.  If  IT  md  Pj  are  used  to  express  and  P^B^  which 

correspond  to  the  (i-l)th  iterative  vector  X*”k,  equation  (7)  can  be 
rewritten  as: 

Xk  -  M^*"1  +  Pk  (8) 

The  following  theorem  is  given  for  diagonal  dominance  cases. 

Theorem  2.5 

If  all  of  the  Jacobian  matrices  in  (1)  are  diagonally  dominant, 
the  Gauss*-Seidel  iterates  of  equation  (8)  converge  to  the  unique 
solution  for  any  initial  point  X° 

Proof 

Assuming  X°  is  an  initial  guess,  by  (7),  we  obtain: 

X1  -  MjX°  +  P± 

X*  -  HjX1  *  P2 

-  HjdljX^Pj)  +  p2 
“  M2M1X°  +  «2P1  +  p2 


k  ^ 


L' 


£ 


K; 


•  •eg 


■  Vk-l ' ' '  ' ' '  "2 VMk-l " '  H3P2 

*Vi-i'"*4V"'*Vk-i4Pt  <»> 

There  exists  s  norm  which  satisfies  | I M I |  -  p(M),  where  p(M)  is  the 
spectral  radius-the  aaxiana  modulus  of  the  eigenvalues  of  the  matrix 
N  [6].  In  ten  of  the  norn,  (9)  will  be: 

lllkll  -  1 I  Vk-1 ' ' '  M2K1  ^Vk-l ' ' '  *2 Wk-l  ” '  M3p2 

♦  Vi-1 — Vs*  •  •  • 4  Vk-i^k  1 1 

By  the  properties  of  the  nora 

IIXk||  <  |lM^Ck_1,“M2M1X°M  +  llMkMk-1‘**M2P1ll 
+  ,,MkMk-l***,l3P2J|  +  I  f*Vhc-l  *  *  ***4P3 1 1 

+  •••  +  HVk-111  +  ,|PkM 

and 

IIXkll  <  ||Mkll  l|Mk„1l|--||M2ll  I  |MXI  I  I  |X°I  I 

♦  11*^11  llMk.1ll***llM2U  llPjl 

♦  m^ii  1 1  Mk__2 1 1  *  ”  1  iMj  1 1  mp2ii 
+  n**kii  iiMk.1!i“*HM4ii  11P3I1 

♦  •••  +  IlMkll  1 IPk„! 1 1  +  llPkll 

suppose  I lMm| J  is  the  aaxiana  norm  aaong  {  IlMjJI  },  and  I |Pm 1 1  <  • 

is  the  aaxiana  nor*  aaong  {  IIPk||  },  the  following  expression  can  be 
obtained: 

Nxkll  <  lluJlMlxjl  ♦  ||pb||  (  1 1^1  |k  +  IIMhII^1 
+  ••••  +  Mill  ) 

Since  the  all  of  the  matrices  JtE^  are  diagonally  doain  nt,  the  spec¬ 
tral  radius  p(M^m^)  <  1  by  Diagonal  Doainance  Theorem  [71. 


Because  p(Ma)  j  l,  so  that  I iM^I I  £  1.  and  it  follows  that 


lia  -  0 
k->«  ^ 

li»  IIXkll  <  ||PBII  (  Hill  -  llM.ll  ) 

k->® 


- 1 


The  Gauss-Seidel  iterative  sequence  {X*}  is  bounded  and  converges  to 
the  unique  solution. 


3  Convergence  Condition  for  Block-Katzenelson-Gauss-Seidel  nethod 


The  motivation  for  adopting  the  block  Gauss-Seidel  iterative 
method  is  to  obtain  convergence  when  the  point  Gauss-Seidel  iterative 
method  does  not  converge,  and  to  improve  the  convergence  rate  (how¬ 
ever,  in  some  cases  the  block  iterative  method  is  slower  than  point 
iterative  method  T9]).  As  long  as  at  least  one  of  the  diagonal  sub- 
matrices  Ajj  contains  a  nonzero  entry  from  above  the  main  diagonal  of 
A,  the  block  Gauss-Seidel  iterative  method  must  have  a  larger  asymp¬ 
totic  rate  of  convergence  than  that  of  point  Gauss-Seidel  iterative 
method  [9] . 

The  Xatzenelson-Gauss-Seidel  method  to  be  discussed  in  this  sec¬ 
tion  is  a  nonlinear  block  Gauss-Seidel  iterative  method  where  the 
blocks  and  subvectors  are  considered  as  basic  units  and  are  solved  by 
Katzenelson’ s  algorihtm. 

In  the  PWL  analysis  of  very  large-scale  systems  and  circuits. 


the  Jacobian  matrices  (J)  of  (1)  can  be  partitioned  into 


where  the  diagonal  blocks  Aii*  . ,n,  which  are  assuaed  to  be 
square  and  nonsingnlar  natrices,  correspond  to  subcircuits,  while  the 
nonxero  eleaents  in  the  off-diagonal  blocks  represent  the  connection 
relationship  between  the  diagonal  blocks  (  snbcircnits). 

The  vectors  1,  X  and  Y  of  (1)  are  partitioned  according  to  the  parti¬ 
tioning  of  Jacobian  natrices;  the  aatrix  equation  (1)  can  then  be 
written  as 


A11  ^2  **•  ^a1 

[,l1 

' 

B1 

A21  ^22  •**  ^a 

X2 

W2 

B2 

F  * 

•  • 

•  • 

* 

• 

+ 

e 

a 

• 

1 

e 

,Anl  ^a2  •••  A*a> 

•V 

•w 

"a 

V 

where  1^,  Wj,  Bj  e  8ni»  i*l....,n  are  the  ith  subvectors 

of  X,  W,  Y  e  Ra  and  ?  n.«  n 

i-1 

The  napping  which  is  the  ith  subset  of  F  can  be  further  expressed 
in  the  fora 

^[(1)  *  ^iiXi+*i  +  Ail^i  “ 


Bi.  j»*i 


(ID 
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Before  we  disease  the  convergence  conditions  of  Xatzenelson- 
Genss-Seidel  method  for  solving  equation  (10)  some  new  definitions 
for  block  Gauss-Seidel  iteration  method  are  needed.  These  defini¬ 
tions  are  created  by  analogy  to  the  point  Ganss-Seidel  case  and  are 
stated  below. 

Definition  3.1 

The  mapping  F:  D  Rn  ->  Rn  is  diagonal  block  isotone  (on  Rn)  if 
for  any  X  6  Rn  the  fractions 

*i,i:  {‘€®1lx+tEie  D}_>  R“l»  #i>i(t)-Fi(X+tEi).  i-1 . m 

are  isotone,  where  E*  is  a  unit  basic  vector  with  ith  snbvector  one 
and  all  others  zero. 

Definition  3.2 

The  mapping  F:  D  Rn  -)  Rn  is  off-diagonal  block  antitone  (on  Rn) 
if  for  any  X  €  Rn  the  fractions 

Pifj:  (t€R1IX+tEj€  D  }->  RnK  j  ( t)-Fi(X+tEj )  i^j .  i.  j-1 . m 

are  antitone,  where  E^  is  a  nnit  basic  vector  with  the  jth  snbvector 
one  and  all  others  zero. 

Definition  3.3 

The  mapping  F:  D  Rn  ->  Rn  is  inverse  isotone  (on  Rn)  if  F(X)  < 

F(Y)  for  any  X,Y  €  Rn  implies  that  X  <y.  Here  X,Y  6  Rn  are  divided 

m 

into  snbvectors  X,  €  R  1.  ial . .  Z  n,-n. 

i-1 

With  these  new  definitions,  we  now  prove  the  following  theorem. 
Theorem  3 .1 


/. 


,  **.  •*.  •*.  **,  v  ■  «  •  .  1 «  '«  *  .  '  *  '  •  *  .  *  .  •  •  *  .  *  • 

v-  v-  -  ~  -V-  -j  -j  -j  v  'j-  v  . 


-  -^-i 


Let  the  PWL  continuous  mapping  F:  D  Rn->Rn  in  (10)  be,  (a)  off- 
diagonal  block  antitone,  (b)  strictly  diagonal  block  isotone,  and  (c) 
inverse  isotone,  then  F  is  defined  as  blook  M-f unction  and  the 
Katzenelson-Ganss-Seidel  method  converges  to  the  unique  eolation. 

Proof : 

The  Katzenelson-Ganss-Seidel  nethod  formula  is  the  following 

. . Xi>-Bi  <12> 

^+1  -  (13) 

where  Ij,  b,  €  Rni,  i*l,...,m  are  the  ith  snbvectors  of  X,  B  €  Rn, 

«  v 

I  n.*n,  and  the  X?  means  the  snbveetor  X.  of  the  k-th  iterate. 
i-1  J  3 

(i)  for  Katzenelson's  algorithm  part 

Here  we  nse  Katzenelson's  algorithm  only  to  solve  the  subcir¬ 
cuits;  that  is,  the  diagonal  blocks,  separately.  Ve  can  consider  that 
equation  (12)  is  the  same  as  eqaation  (11).  Let  then  (11)  can 

be  rewritten  as: 

Fi<xi>  ■  xiV»i  -  ci 

ci  * 

When  the  Ganss-Seidel  iterative  method  is  applied  to  the  matrix  equa¬ 
tion  such  that  X*»(X*,. . , ,X*)^  and  X*+*,  j“l,...#i-l  are  given  from 
previous  Gaass-Seidel  iterate  and  the  first  i-1  step  in  this  iterste, 

so  that  ia  (li)  is  considered  constant.  Then,  F^  is  a  piecewise- 
linear  mapping  corresponding  to  subcircnit  i,  and  Katzenelson's  algo¬ 
rithm  is  applied  to  solve  the  subcircnit. 
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There  are  already  a  number  of  paper*  dealing  with  the  conver¬ 
gence  properties  of  Eatzenelson’ a  algorithm.  We  jnst  mention  some 
conclusions  here.  Katzenelson* s  algorithm  was  originally  applied  to 
the  nonlinear  network  consisting  of  resistors  with  continuous, 
piecewise-linear,  strictly  isotonic  characteristics  [1].  Later,  it 
was  proved  that  if  the  determinants  of  in  (1)  had  the  same  sign 
in  all  the  regions,  then  Eatzenelson' s  algorithm  was  guaranteed  to 
converge  to  a  solution.  In  this  paper  we  only  deal  with  the  circuits 
which  consist  of  the  nonlinear  resistors  with  continuous  piecewise- 
linear  strictly  isotonic  characteristics  so  that  the  Eatzenelson' s 
algorithm  is  guaranteed  to  converge  to  a  unique  solution. 

(ii)  for  Gauss-Seidel  part: 


Suppose  the  starting  points  are  X°,  T°  €  Rn,  and 


0\T 


v°  z  (x°  T° 

A 

Y°  >  (Y°  Y°  Y°)^ 

and  V*  is  the  solution, 

V*  *  (V*  V*  vV 

where  X^.Y^  and  V*,  i=l,...,m  are  subvectors  corresponding  to  sub¬ 
circuit  i.  Since  F  is  block  isotone,  the  following  relations  are 
satisfied 

F(Y°)  >  p(V*)  >  F(X°).  means  that.  Y°  >  V*  >  X° 


By  induction,  for  some  X  >  0,  i  >  1, 

X°<Xk<V*<Yk<Y°,  F(Xk)  <  B  £  F(Yk)  (14) 

<  xf 1  <  v;  <  rf 1  <  yJ.  j-i . i-1  (15) 


For  k-0,i»l,  (14)  hold*,  and  j-0  in  (15),  *o  that  (15)  is  empty.  Snp- 
poae  (14)  and  (15)  hold  for  k  and  i-1. 

By  the  off-diagonal  block  antitone  property,  it  follows  that 


P(S) 

“  Fi(xJ+1 

$  *  •  e  0 

S'Xi+l'*‘ 

.) 

(15) 

Q(S) 

m  F  (Yk+^ 

i|Xl 

$  •  •  •  9 

S,*i+1‘*  * 

.) 

(17) 

then. 

Q(S)  < 

P(S) 

S  € 

8ni 

(18) 

and 

P(Xk) 

<  F.(xk) 

iBi 

i  F.(Yk) 

<  Q(Yi> 

P(VJ) 

>  B.  -  F 

i(vj. 

y* 

..V*)  >  Q(vJ) 

(19) 

Let  P(Xk)  «  Bj  ■  Q(Yk) ,  and  by  the  diagonal  block  isotone  and  (19), 
we  can  obtain 


and 


^  x  rj 


For  Ganss-Seidel  method 


Yk+1  -  Yk 

i  i 


Xk+1  -  Xk 

4i  4i 


then. 


xk  <  xk+l  <  y*  <  yk+l  <  yk 


holds  for  i*l,2,...,a,  and  hence, 

xk  <  yk+1  <  y*  <  yk+1  <  yK^ 

Then,  we  obtain, 

F.(Yk+1)  >  F.d^1 . Yk+1.Y k....) 
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Similarly* 

Fi(Xk+1)  <  B£ 

Hence,  this  completes  the  induction  of  (14)  and  (15),  and  the 
sequences  (X)  and  {Y}  obtained  from  Ganss-Seidel  iterate  have  the 
limits 

lim  Xk+1  -  X* 
k->“ 

and 

lim  Y^1  -  Y* 
k->“ 

and  the  following  formulas: 

F(X*)  <  B,  and  F(Y*)  >  B 

By  the  continuous  Mapping  ,  it  follows  that  F(X*)»F(Y*)“B, 

.  _•  •  _• 

and  X  »Y  «V  . 

Remark  3.1 

In  the  above  proof,  diagonal  block  isotone  means  the  mapping 
is  isotone  in  the  diagonal  block  i  with  respect  to  subvector  X^. 
Actualy,  it  is  not  necessary  that  the  mapping  F^  be  isotone  with 
respect  to  the  entire  subvector  X^.  For  instance,  assume  the  mapping 
F^  in  (11)  can  be  partitioned  into  two  components  F-^and  F^: 


il^r^r  + 

Ar.Xs 

Z  A  ,x. 

j-lrJ  J 

*  Wr  -  Br 

(20) 

i2*AsrXr  + 

A..X. 

£  AsjXj 
J-l  3  J 

+  W.  “  B. 

(21) 

j-l,....m,  j4i 


where  Xy,  Xf  and  Br,  Bs  are  the  components  of  Xj  and  B^: 

Xi  “  <Xr.X,)T.  Bi  -  (Br.B#)T* 

■“d  is  partitioned  into: 


V.V.Vi 


and  Xi#  Bj  e  Rni.  Xr,  Br  €  Rnr,  X#,  B#  €  Rns.  nr+n#-ni. 

Assume  also  that  X#  related  to  the  other  subvectors  through  the 
off-diagonal  faloeks  and  while  Xy  it  not  related  to  the  other  subvec¬ 
tors,  that  means  A^,  j«l,...,m,  jj*i,  are  zero,  and  (20)  becomes: 

ArA  *  Ar.*.  +  »r  *  Br 
so  that  Xf  can  be  expressed  as: 

\  ‘r,1.  ♦  C'W  <“> 

Substituting  (22)  into  (21),  we  obtain: 


where 


K-hi1.  *  ♦*,  -  b;.  i*i 


A  ii*  Ass  ~  ^sr  Arr  Ars 


B'i  -  B,  -*rr  <Br  ‘  V 

We  can  now  nse  (23)  instead  of  (11),  during  the  Ganss-Seidel  itera¬ 
tion.  Then,  it  is  sufficient  to  have  isotone  with  respect  to  sub- 
vector  X#  instead  of  X^  for  convergence. 

Remark  3.2 

It  is  not  necessary  that  the  F  be  inverse  isotone  with  respect 
to  the  entire  vector  X  either.  In  order  to  explain  this,  we  still  use 
the  above  equation  (20),  (21),  (22),  and  (23).  By  eliminating  Xf 
that  is  not  related  to  the  other  snbvectors  through  the  off-diagonal 
blocks,  matrix  A^^  becomes  A'^j  and  snbvector  Xj  becomes  X'j“Xs  and 
becomes  F'^.  Similarly,  the  other  subvectors  Xj ,  j«l,...,m. 


•  ■  • 


be cone  X’j,  j-l 
a  new  mapping  F 


.  ,m,  and  the  diagonal  blocks  Ajj  become  A'jj 

t 

is  obtained  which  consists  of  F.,  ixl,...,m. 

t 


Then 
It  is 


sufficient  that  mapping  F  be  inverse  isotone  with  respect  to  subvec¬ 
tors  X'j  i-1 , . . . ,m  for  convergence.  Using  X'^  and  A'jj  instead  of  Xj 
and  A^^(  namely,  using  the  mapping  F'  instead  of  the  mapping  F, 
Theorem  3.1  still  holds. 

Remark  3.3 


It  is  possible  that  the  convergence  of  the  Katzenelson-Gauss- 
Seidel  method  can  be  speeded  up  by  using  SOR  iteration  instead  of 
standard  Gauss-Seidel  iteration,  that  is,  using 

Xi+1  "  Xi  +  “<Xi  "  XX)»  0  <  «  <  2.  (24) 
instead  of  equation  (4). 

Another  factor  which  affects  the  convergence  rate  is  the  abso¬ 
lute  values  of  the  off-diagonal  elements.  The  smaller  the  absolute 
values  of  the  upper-triangular  elements  in  the  whole  matrix, the  fas¬ 
ter  the  convergence  will  be.  This  proof  can  be  found  in  [7]. 

So  far,  we  have  only  discussed  the  properties  of  mapping  F.  We 
show  now  what  kind  of  matrices  will  ensure  that  F  in  (10)  is  block 
M-function. 

Definition  3.4 


A  matrix  A  is  block  M-matrix  if  the  matrix  A  can  be  partitioned 
into  blocks,  such  that  (a)  all  of  the  elements  in  off-diagonal  blocks 
are  negative  or  zero,  (b)  each  diagonal  block  is  diagonally  dominant 
with  each  diagonal  element  positive  (c)  each  diagonal  block  can  be 


tranformed  into  a  class  M-matrix  by  means  of  the  method  in  Remark  3.1 


and  3.2 


Theorem  3.2 

If  the  natrices  {J}  in  (1)  or  in  (6)  are  block  M-matrices,  then 
the  Katzenelson-Gauss-Seidel  method  it  convergent. 

Proof: 

Clearly,  fron  the  condition  (a)  and  (b)  in  definition  3.4,  it 
follows  that  F  is  off-diagoanl  block  antitone  and  diagonal  block  iso¬ 
tone.  According  to  condition  (c),  each  diagonal  block  becomes  M— 
matrix  after  some  elements  are  eliminated.  The  inverse  of  an  M-matrix 
is  a  positive  matrix  which  ensures  the  diagonal  blook  to  be  inverse 
isotone  with  respect  to  those  variables  which  have  relationships  with 
the  other  diagonal  blocks.  Since  the  diagonal  elements  are  dominant 
in  the  entire  matrix  and  the  off-diagoanl  elements  are  negative,  the 
entire  matrix  is  an  H-matrix  after  each  diagonal  block  is 
transformed  into  an  M-matrix  by  means  of  eliminating  the  variables 
that  have  no  relationship  with  the  other  blocks.  Then  with  the 
Remark  3.2,  the  mapping  F'  is  inverse  isotone  with  respect  to  the  new 
snbvectors  and  F’  is  block  M-  function  and  the  Katzenel son-Gauss- 
Seidel  method  is  convergent. 

The  other  case  is  when  the  entire  matrix  and  its  diagonal  blocks 
are  not  diagonally  dominant,  then  by  using  the  transformation  formu¬ 
las  in  Remarks  3.1  and  3.2,  the  new  reduced  matrix  and  its  diagonal 
blocks  become  diagonally  dominant  and  the  block  Gauss-Seidel  method 
becomes  convergent. 

The  following  are  examples  where  the  block  Gauss-Seidel  method 
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converges,  while  the  point  Gauss-Seidel  does  not. 
Example  1 

Consider  the  matrix  equation: 


'3  7  0  0 

8  4  0  -3 

X 

x1 

2 

’  5 
4 

0  0  5  -3 

X3 

8 

.0  -2  8  6, 

h3J 

l10J 

If  the  point  Ganss-Seidel  method  is  nsed  to  solve  this  equation,  the 
following  formula  is  obtained: 


fjll 

k+l 

f3 

0 

0 

0 

-y 

'0 

-7 

0 

0 

X2 

8 

4 

0 

0 

.  f 

0 

0 

0 

3 

X3 

0 

0 

5 

0 

0 

0 

0 

3 

K 

lo 

-2 

8 

6 

V 

0 

0 

0 

0 

in  which  case  p( (D-L)-1U)=4.6  and  the  iterates  will  not  converge.  If 
the  block  Gauss-Seidel  method  is  adopted  to  solve  this  equation  and 
the  4x4  matrix  is  partitioned  into  2x2  blocks,  the  following  formula 
is  obtained: 


,k+l 


i  X, 
1*4 


3  7  0  <h 

8  4  0  0 

0  0  5  -3 

U>  -2  8  6 


rl 


KO  0  0  0 

0  0  0  3 

0  0  0  0 

10  0  0  01 


rx. 


sX4 


The  spectral  radius  of  the  bolck  Gauss-Seidel  matrix  is  0.341,  which 
is  less  than  1,  so  that  the  block  Gauss-Seidel  method  is  convergent. 
The  new  reduced  matrix  generated  from  the  transformation  described 
above  is 


-14.66  -3  ' 

.-2  10.8  , 


which  is  diagonally  dominant. 


Example  2 


Given  the  matrix: 

f3  -7  0  0 

8  4  0  -3 

0  0  5  -3 

>0  -2  8  6 


It  can  be  transformed  into: 


(  22. 

I"2 


66  -3 

10. 


by  using  block  Gauss-Seidel  method.  The  reduced  matrix  is  a  class 
M-matrix  and  a  matrix  equation  consists  of  this  matrix  can  be  solved 
by  block  Gauss-Seidel  iterative  method. 

In  general  the  Katzenelson-Gauss-Seidel  method  can  be  used  on  a 
wider  class  of  problems  than  having  the  reduced  matrix  digonally  dom¬ 
inant  or  class  M  matrix.  For  simplicity,  we  only  discuss  the  cases 
in  which  matrices  are  partitioned  into  4  blocks.  For  instance,  A^, 
Aj,  A3,  A^  are  blocks  of 

the  main  matrix  A: 

.  _  fAl  M 


A4i 

Aj^A^  are  diagonal  blocks,  and  Aj,  A3  are  off-diagonal  blocks. 
Consider  a  linear  equation  consisting  of  the  matrix  A  and  the  vectors 

X'i'*X1,X2*T'  *8  follo,r8: 


At  *2 

B1 

.A3  A4 

.*2 

“ 

•®2 

Applying  block  Gan* s-Se idol  method,  ve  obtain  the  foil owing: 
k+1 

_  f*i 


*i 

X 


*1  Of*  f°  -*2]  _  “f  fBl 

*J  io  oJ  lljl+[*3  *J  Uj 


rk+1 


MXk  +  P 


or,  X 

where  H  and  P  are  respectively: 


(25) 

(26) 


Theorem  3.3 

If  the  maximum  modulus  of  the  eigenvalues  of  matrix  N  in  (26)  is 
smaller  than  1,  the  block  Gauss-Seidel  iterative  process  converges  to 
the  unique  solution. 

Proof: 

Assuming  X°  is  initial  guess,  by  (26),  we  obtain: 

X1  =  MX°  +  P 

X2  -  MX1  +  P  -  M(MX°  +  P)  +  P 

-  P  +  MP  +  M2X°, 

•  •  •  9 

Xk+1=  MXk  +  P 

=  P+  MP  +  M2P  +  ...,+  Mk-1P  +  MkX° 

-  P(I  +  M  +  M2  +  ...,+  Mk_1)  +  MkX° 

By  Newman  Lemma  [6],  if  p(M)  <  1,  then 

lim  T  Mk  *=(  I-M) -1 
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lia  ■  0 
k->« 

such  that 

lia  x^1  -  P(I-M)"1 
k->« 

Jr  1  converges  to  the  unique  solution. 


The  matrix  M  is  the  product  of 

f*i  of*  ro  -AJ 


and  beooaes 


[At  °]  [°  "*2 
1*3  V  l*0  0  ■ 

M-f°  ^ 

Lo  a"1a3a^1a2. 


This  is  equivalent  to  requiring  that  the  aaziana  modulus  of  the 
eigenvalues  of  A  ^AgA^Aj  is  less  than  1.  This  is  the  basic  conver¬ 
gence  condition  for  block  Gauss-Seidel  method  in  the  linear  case.  In 
this  case,  the  matrix  A  need  not  necessarily  be  of  class  K-matrix  or 
diagonally  dominant. 


In  piecevise-1 inear  cases,  the  elements  of  A^,  Aj,  Ag,  and  A^ 
will  change  their  values  according  to  the  values  of  the  variables.  If 
the  block  Katzenelson-Gauss-Seidel  method  is  used,  equation  (26) 
should  be  rewritten  as: 

Xk+1  -  \+1Xk  +  Pk+1  (27) 
where  end  Pj.  are  variable  depending  on  the  piecewise  linear  region 
k  which  the  vector  X*  is  located  in. 
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Theorem  3.4 

A  piecewise-linear  mapping  P  in  (1),  or  (10)  can  be  split  and 
expressed  as  in  (27)  when  block  Ganss-Seidel  iterative  method  is 
applied;  if  the  spectral  radii  pd^)  tre  less  than  1,  k«l,...,r,  then 
Katzenelson-Ganss-Seidel  method  converges. 

The  proof  is  similar  to  that  of  Theorem  2.5  except  that  the  matrices 


'H^'s  are  different.  The  is  the  block  Ganss-Seidel  matrix  here* 
while  the  is  theorem  2.5  is  the  point  Ganss-Seidel  matrix. 


4  Examples 


Example  1 

The  first  example  is  a  bipolar  circuit  which  is  shown  in  Fig. 5, 
whose  snbcircuits  are  TIL  NAND  gates  as  the  one  shown  in  Fig.  1.  The 
bipolar  transistors  are  expressed  by  Ebers-Koll  model  shown  in  Fig. 2, 
and  the  PWL  function  in  Fig. 3 (a)  represents  approximately  the 
exponential  function  for  diode  characteristics  in  Ebers-Moll  model. 
Then  the  PWL  model  which  replaces  each  npn  bipolar  transistor  in  PWL 
analysis  is  obtained  in  Fig. 3(b).  The  PWL  characteristics  is 
strictly  monotonic  and  expressed  by  seven  segments,  and  each  segment 
is  described  by  four  parameters  (start  breakpoint,  slope,  intercept 
and  endbreak  point).  The  parameters  of  the  diode  characteristics  are 
given  in  table  1. 

Table  1 


segment 

end  breakpoint 

slope 

intercept 

1 

-0.1 

0.1755e-10 

-0.4580d-7 

2 

0.1 

0.4581 e-6 

0.0e+0 

3 

0.3 

0.5127e-3 

-0.5122e-4 

4 

0.5 

0.1124e+l 

-0.3371e+0 

5 

0.7 

0.2462e+4 

-0.1231e+4 

6 

0.9 

0.5396e+7 

-0.3777e+7 

7 

m 

0.1182e+ll 

-0.1064e+ll 

where  the  slopes  and  intercepts  are  taken  respectively  as  conduc¬ 
tances  and  independent  current  sources  in  the  equivalent  circuits. 
The  parameters  a^o.99,  <*r«0.5  for  the  Ebers-Moll  model  in  Fig. 2  ere 
considered  constant. 

The  circuit  in  Fig.S  is  partitioned  into  4  suboireuits  as  shown 
in  Fig. 6(a),  where  small  resistors  are  inserted  in  such  a  way  that 
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each  subcircuit  input  and  ontpnt  nodes  will  not  be  connected  to  other 
snbcircnits  directly:  rather  they  are  connected  to  a  set  of  new 
nodes:  nl  ~  n5.  Since  the  resistors  inserted  are  the  sane,  the  rela¬ 
tions  of  new  node  voltages  to  the  snbcircnit  input  and  ontpnt  node 
voltages  are 

Vnl“( v7+v>  9+^8) /3 ' 

Vn2“(v8+vV/2» 

▼^-(vj+Vj)/!,  (28) 

Tn5"(W/2‘ 

where  ~  y^  belong  to  first  snbcircnit,  v'^  ~  v'g  belong  to  the 
second  snbcircnit,  v^  ~  y^  belong  to  the  third  snbcircnit,  and  v^  ~ 
belong  to  the  fourth  snbcircnit.  The  small  resistors  are  chosen 
to  be  0.1  ~  0.001  ohm  so  that  the  conductances  are  10  ~  1000  mhos. 

When  the  entire  circuit  reaches  the  equilibrium  point,  v^wy^wv' ^«*Vg , 

Tn3*VV  Tn4"V**V  vn5*W 

The  matrices  of  the  partitioned  networks  with  the  inserted  resistors 
are  formulated  by  the  modified  nodal  approach  and  their  structure  is 
shown  in  Pig.  7. 

When  the  Katxenelson-Gauss-Seidel  method  is  applied  to  solve  the 
matrix  equation  of  this  bipolar  transistor  circuit,  the  iterative 
procedure  is  as  follow: 

(a)  nl  ~  n5  are  taken  as  tearing  nodes; 

(b)  initial  guesses  for  v^  ~  y^  are  given; 
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(c)  each  subcircuit  it  separated  such  that  each 
subcircuit  appears  to  be  drives  by  independent 
voltage  sources  v^  -  as  in  Fig.  8; 

(d)  Katzenelson's  algor it ha  is  applied  to  solve  each 
snbeircnit  at  a  set  of  gnesses  of  v^  ~  y^* 

(e)  nse  the  solution  to  the  subcircuit  to  update 
the  set  of  v^  -  y^  by  (28),  and  go  to  solve 
next  snbcircuit  by  Katxenelson' s  algorithm; 

(f)  after  solving  all  of  the  suboircuits,  obtain 
new  set  of  v^  -  Vji5  in  tens  of  (28); 

(g)  check  the  new  set  of  v^  ~  with  old  set  of  guesses, 
if  their  differences  are  less  than  a  given  tolerance, 
stop  iterating,  otherwise,  take  the  new  set  of  v^  -  y^ 
as  guesses  and  go  to  the  step  (d)  and  iterate  again. 

As  was  discussed  above,  if  the  new  reduced  diagonal  blocks, 
which  are  generated  froa  the  diagonal  blocks  by  eliainating  internal 
node  voltage  variables  in  the  corresponding  suboircuits,  are  diago¬ 
nally  doainant  or  class  M-aatrices,  then  the  Katzenelson-Gauss-Seidel 
method  converges.  The  following  are  some  examples  whioh  show  how  to 
generate  the  new  blocks. 


✓4 
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The  nodal  matrix  of  the  snbcircnit  1  in  a  given  region  is: 


0000 

-24.87 

4975 

-2462 

0000 

0000 

-2462 

0000 

0000 


0000 

-.4e-6 

-,4e-6 

2488 

-2462 

-24.87 

0000 

0000 

0000 


0000 
0000 
0000  - 
-24.87 
-2462 
4949  - 
-2462 
0000 
0000 


-24.8  - 
-2462  - 
0000 
0000 
0000 
0000 
0000 
le+10 
0000 


.011  ’ 
1.12 
0000 
0000 
0000 
0000 
0000 
0000 
le+10, 


In  the  snbcirenits  (or  blocks).  v^~v^  ut  not  related  to  the  other 
snbcircnit  (or, blocks),  they  are  only  internal  node  voltages  which  do 
not  appear  in  the  Ganss-Seidel  iterative  process.  Using  the  method  in 
Remark  3.1,  we  can  obtain: 


4950 

4924 

0000 

0000 

0000 

0000 

-7385 

9873 

-24.87 

-.4e-6 

0000 

0000 

0000 

-2480 

4975 

-,4e-6 

0000 

0000 

0000 

2462 

-2462 

2488 

-2462 

-24.87 

0000 

0000 

0000 

-2462 

4934 

-2462 

0000 

0000 

0000 

-24.87 

-2462 

4949 

0000 

-24.8 

-.011' 

0000 

-2462 

-1.12 

-2462 

0000 

0000 

0000 

0000 

0000 

0000 

0000 

0000 

-2462 

0000 

0000 

0000 

0000 

-2462 

0000 

0000 

-2462' 

-24.8 

-2462 

0000 

0000 

0000 

0000 

.  2461 

-2462 

0000 

0000 

0000 

0000, 

7486 

0000 

0000 

* 

0000 

le+10 

0000 

.  0000 

0000 

le+10 

4 

and  the  new  block  is  A' ii“A,f"AsrAr,  Art,  which  is 


f  ,261e+4  -,103e+2  -.471e-2  ’ 

-.122e+2  .100e+ll  -.113e+l 

l-.667e-l  -.119e+l  .lOOe+11- 


This  is  *  class  M-matrix  as  well  as  diagonally  don inant. 


The  nodal  matrix  of  snbcircnit  3  in  a  given  region  is 

4924  4924  0000  0000  0000  0000  0000  0000  -.17e-12 

-9847  9848  -.46e-8  -.17e-10  0000  0000  0000  0000  -0.17e-10 

17e-10  0000  0000  -.17e-10 

2487  -2462  -24.87  0000 

2462  4934  -2462  0000 

24.87  -2462  2486  -.51e-3 

0000  0000  -.5 le-3  100.0 


which  can  be  transformed  into: 


.100e+3  .24e-21  .15e~10 

,666e-16  .100e+3  -.853e-l 

,645e-15  -,136e-10  .100e+3 


This  new  block  is  diagonally  doaiinant. 


Ve  can  also  transform  the  other  blooks  into  new  blocks  which  are 
diagonally  dominant  and  sometimes  are  also  M-matrices,  so  that  the 
Jacobian  matrices  of  the  entire  circuit  are  M-matrix  or  diagonally 
dominant  also.  The  spectral  radii  of  the  block  Ganss-Seidel  matrices 
are  less  than  1,  and  the  Katzenelson-Gauss-Seidel  method  is 
guaranteed  to  be  convergent. 

Note  that  in  the  original  circuit  in  Fig. 5,  there  are  no  resis¬ 
tors  connecting  the  subcircuits  together,  and  v'-j^Vg, 


v-»*Vft.  Since  small  resistors  are  inserted  between  the  subcircuits. 


the  accuracy  of  node  voltagea  is  affected.  We  take  the  absolute 
differences  between  v^,  v*g  and  Vg,  and  the  difference  of  v'^  and  Vg, 
and  the  difference  of  v^  and  Vp  as  the  accuracy  measurement,  for 
instance,  v^-3.99,  y'p=3.89,  the  accuracy  is  0.1.  The  following 
table  discribes  the  relation  between  the  conductance,  the  iterative 
number  and  the  accuracy. 

Table  2 


conductance 

iterative  number 

accuracy 

10 

13 

0.1 

100 

33 

0.01 

1000 

165 

0.005 

■If  fewer  resistors  are  inserted  as  shown  in  Fig. 6  (b),  faster 
convergence  rate  is  obtained,  since  the  sabcircait  outputs  become  the 
'independent  voltage  sources',  which  are  now  applied  to  directly 
drive  the  next  subcircuits. 

Example  2 

In  this  example  which  is  an  MOS  circuit  there  is  no  need  to 
insert  small  resistors  in  order  to  apply  the  Katzenelson-Gauss-Seidel 
method.  Instead  of  attempting  to  obtain  a  bordered-block-diagonal 
matrix  structure,  as  is  dine  in  Example  1,  the  circuit  is  decomposed 
into  dc-connected  subcircuits  and  then  analysis  sequence  is  performed 
to  order  the  equations  into  nearly  power-block-diagonal  form  [12] 
The  example  circuit  shown  in  Fig. 9  is  an  NOS  register  circuit.  In 
this  circuit  and  in  the  following  example  circuit,  the  NMOS  transis¬ 


tor  will  be  used  to  illustrate  the  method.  The  transistor  model  is 
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not  discussed  here  in  detail;  the  mathematical  equations  for  the 
channel  current  in  the  NMOS  transistor  can  he  simply  expressed  as 
til]: 

Ids“K<Vgs-Vt>2  "  K{Vgd"Vt>2  (29) 

where  Id#  is  the  channel  current  from  Source  to  Drain,  Vgs  is  the 

voltage  across  Gate  and  Source,  Vg<J  is  the  voltage  across  Gate  and 
Drain,  is  the  threshold  voltage.  This  yields  the  Ebers-Moll 
'like'  model  of  Fig. 10(a).  Both  forward  and  reverse  a  must  equal  to 
unity  in  order  for  the  gate  current  Ig»0.  as  it  must  be  for  an  NOS 
transistor.  Fig. 10(b)  shows  the  HfL  characteristics  which  is  the 
approximation  of  the  quadratic  characteristics  for  a  driver  NOS 
transistor.  The  FffL  characteristics  for  a  load  NOS  transistor  is  the 
same  shape  as  that  of  the  driver  except  that  the  origin  point  and  the 
vertical  axis  are  shifted.  For  digital  applications  it  was  fonnd 
that  only  three  segments  are  sufficient  to  describe  the  characteris¬ 
tics.  The  parameters  for  PWL  characteristics  are  as  follows: 

Table  3  for  driver  transistor 


segment 

breakpoint 

slope 

intercept 

1 

1.0 

1.0e-12 

-1.0e-12 

2 

2.25 

18.0e-5 

-18.0e-5 

3 

m 

82.e-6 

-154. e-6 

Table  4  for  load  transistor 


segment  breakpoint  slope 

1  -2.0  1.0e-12 

2  -1.0  6.25e-6 

3  «  18.75e-6 


intercept 

2.0e-12 

12.5e-6 

25.0e-6 


The  entire  circuit  is  divided  into  3  snbcircnits,  as  in  Pig. 11, 
and  v j  are  the  input,  Vj  ~  vg  are  to  be  solved.  Ve  first  solve 
the  inverter  in  Fig. 11(a),  and  obtain  vg  which  is  the  solution,  since 
v^  is  affected  by  Vj  only.  For  the  subcircnit  in  Fig. 11(b),  v^,  vj 
and  Vg  are  given,  Vg  is  not  known  and  has  to  be  guessed;  then  v^  and 
v3  are  solved  in  tens  of  Katzenelson's  algorithm.  After  v^  is 
obtained,  it  becomes  the  input  to  the  inverter  in  Fig. 11(c)  which  is 
also  solved  by  Katzenelson' s  algorithm.  The  new  v^  i,  compared  with 
the  guess;  if  they  are  the  same,  the  solution  is  obtained;  otherwise, 
go  to  the  subcircuit  (b),  and  start  the  next  iteration,  until  conver¬ 
gence  is  obtained.  In  this  example,  the  Jacobian  matrices  of  the 
entire  circuit  are  class  K-matrix  in  most  of  the  regions,  but  not 
always.  For  instance,  when  the  initial  guesses  are: 


-  0.0, 

▼2  *  °*o. 

▼3  -  0.305. 

▼4  ”  0.38e-8, 
v5  »  5.0, 

▼g  -  0.305, 

the  matrix  is: 


'  0.20e-ll 
0.0000 
0.0000 
l 0.0000 


0.00000 

0.82e-4 

-0.10e-ll 

0.00000 


0.0000 

-0.10e-ll 

0.20e-ll 

0.00000 


0.0000  ■ 
0.0000 
0.0000 
0.18e-4 


When  they  reach  the  solution: 
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-  0.0, 

▼2  *  0.0, 

▼3  -  0.305, 

▼4  -  0.305, 

▼5  -  5.0, 

▼g  *  5.0, 

the  Jacobian  matrix  ia: 

'  0.18e-4  0.00000  0.0000  0.0000  ' 

0.0000  0.16e-3  -0.82e-4  0.0000 

0.0000  -0.82e-4  0.82e-4  0.0000 

L 0.0000  0.00000  0.00000  0.18e-4' 

During  the  Gauss-Seidel  iteration,  the  matrix  elements  change  their 

values  but  they  keep  the  same  signs  and  diagoanal  dominance  is 

preserved,  so  that  the  Jacobian  matrices  are  class  K-matrices. 


The  Jacobian  matrices  are  not  class  M*-matrix  in  some  PYL  regions 
during  Katzenelson’ s  process,  when  other  initial  guesses  are  used; 
for  instance,  vhen  the  subcircuit  (b)  is  been  solving  the  node  vol¬ 
tages  are: 
v2  -  0.0, 

▼2  -  5.0, 

▼3  -  3.0, 

▼4  *  0.36e-7, 

▼5  -  5.0. 

▼g  =  0.305, 

the  Jacobian  matrix  is: 

’  0.82e-4  0.00000  0.0000  0.0000  ' 

0.0000  0.18e-4  -0.10e-ll  0.64e-4 

0.0000  -0.10e-ll  0.82e-4  0.0000 

1 0.0000  0.00000  0.00000  0.18e-4' 

which  is  not  class  M-matrix 

After  the  Katzenelson’ s  process  of  subcircuit  (b)  reaches  the  solu¬ 
tion:  Vj  »  0.305,  v4  *  0.37e-8,  with  the  other  voltages  the  same,  the 
Jacobian  matrix  becomes  class  M-matrix  again: 
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' 0.82e-4  0.00000  0.0000  0.0000 

0.0000  0.82e-4  -0.10e-ll  0.0000 

0.0000  -0.10e-ll  0.82e-4  0.0000 

U.0000  0.00000  0.00000  0.18e-4< 

From  the  viewpoint  of  the  Gauss-Seidel  iterations,  the  Jacobian 
matrices  are  class  M-matrix  in  the  region  where  the  solution  exists 
and  the  Katzenelson-Gauss-Seidel  method  converges  although  the  Jaco¬ 
bian  matrix  may  not  be  class  M-matrix  in  some  regions. 

Example  3 

The  third  example  is  a  large-scale  MDS  circuit-a  programmable 
logic  array  (PLA),  which  is  used  to  implement  a  traffic  light  con¬ 
troller  [13].  The  circuit,  which  is  shown  in  Fig.  12,  is  partitioned 
into  42  subcircuits.  Some  feedback  node  voltages  are  assumed  to  be 
given  by  previous  iterates,  according  to  the  Gauss-Seidel  method.  The 
subcircuits  are  solved  by  Katzenelson  algorithm  in  a  sequence  which 
is  obtained  by  topological  properties.  After  three  iterations,  the 
solution  is  reached.  The  total  computational  time  is  0.683  second. 
This  example  shows  that  the  Katzenelson-Gauss-Seidel  method  is  suit¬ 
able  for  large-scale  circuits. 
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Appendix 


Appendix  1 

The  following  two  proofs  can  be  found  in  [6],  which  are  repeated  here 
for  easy  reference. 


Proof  of  Theorem  2.2 


Supposing  starting  points  are  X°,  Y°  €  Rn, 

t°„ / _o  o  Ot  T  *o  /  o  o 

^  *  *•  **y^» •  •  •  »7i # • 

and  mapping  F  can  he  expressed 
MX)  -  b. 


y°)T 


*£(X)  -  b£ 


As  induction  hypothesis. for  some  X  >  o,  i  i.  1, 

X°  <  Xk  <V*  <  Yk  <  y°  F(X)k  <  B  i  F(T)k  (30) 

*k<[x J+1  <y J*1  <y j  •  J  -1 . 1-1 .  (31) 


For  k*0,i*l,  (30)  holds,  and  (31)  is  empty,  suppose  (30)  and  (31) 
hold  for  k  and  i— 1.  By  the  off-diagonal  antitone  property,  it  fol¬ 


lows  that  the  functions 

P(»Wi(x^1.....xi_k+1,s,x1+k . xj)  (32) 

Q(s)-fi(yX+1,...,yi_J+1,s,yi+J,...,yJ)  (33) 

satisfy 

Q(s)  <  P(s)  ,  and  s  e  R1  (34) 
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Q(*k)  <  P(xk)  <  £l(Xk)  <  bi  1  fi(Yk)  <  Q(yJ)  <  P (y\) 

Let  P(xk)  ■  bj  ■  Q(yk),  by  the  strictly  diagonally  isotone  property. 


we  obtsin 


xk  /  -k  ,  -k  .  k 
xi  i  *i  i  y£  i  yj. 


where  the  relation  x*  <  yk  is  the  consequence  of  (34) 

For  Gauss-Seidel  method  :  y^+*  yi 

xf  1  -  3 


then,  xk  <  xk+1  <  yk+k  1  yk  holds  for  i-l,2,...,n  and  hence, 
Xk  <  Xk+1  <  Y^1  <  YK. 


Then,  we  obtain. 


f i(Yk+1)  >  fi(yk+1 . 7i+1.yi+k»...)“bi. 


FCY^1)  >  B 


Similarly, 


F(Xk+1)  <  b 

This  completes  the  indnction  and  hence  the  proof  of  (30)  and  (31). 

Clearly,  now  the  limits 

limXk-X*  <  Y*»li«Yk 
k->o»  k->« 

exist  and  by  the  definition  and  the  continnons  mapping  ,  it  follows 
that  F(X) *=  F(Y)*»  B,  and  X*  -  Y* 


Proof  of  Theorem  2.3 


For  given  X°,B  €  Rn.  define 

D°  -  F“1(«extf1(X0),B1],...,**x[fn(X°).Bn]) 


V°  -  F"1(«intf1(X°).B1] . ain[fa(X°)  »BnJ ) 

By  the  inverse  isotonic ity 

F(U°)  >  B  l  F(V°) .  means  that  U°  >  X°  >  V°  end  D°  >  X*  >  V° 
Let  {Uk},  {Vk}  and  {Xk}  denote  the  Gauss-Seidel  iterative  seqnenoes 
starting  fro*  D°,  V°  and  X°  respectively.  By  the  property  of  strictly 
diagonally  isotone  and  continuity  of  F.  the  solutions  tL,  and  x. 

of  equations 


f  AtH?1.. 

...u.  k+k  Sk  n.  k. ... 

,  uk 

i-1  ' 

'  i+1 - 

*  n 

fl(vk+1.., 

k+1  -k  „  k 

,vk 

i'xl 

i-1 

i+1 - 

*  n 

fl(xk+1... 

„  k+1  -k  „  k 

,xk 

ivxl  ' 

i-1  ' 

i+1 - 

n 

i-1... 

» • I&I  k«0 #1 

M  • 

exist  and  are  unique  and  therefore  the  Gauss-Seidel  sequences  are 
well  defined.  Considering  U°,  V°  as  X°  and  7°  in  theorem  2.2,  we  can 
obtain 


y°  <  yfcfl  <  yk  <  liByk  <  liaDk  ^  Dk+1  <  Uk  <  U° 

k->-  k->« 


F(Vk)  <  B  i  FiV*) 

For  {X},  suppose  that  for  sone  k  >  o,  i  2  1 

Vk  <  Xk  <  Uk.  vk+1  <  ,k+1  <  uk+1.  j-1  . i-1 
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Append ix  2 

The  PLA  eirenit  is  solved  by  neans  of  the  PREMOS  program  [12] 
vhich  is  modified  to  inolnde  the  PWL  analysis  part.  The  following  is 
the  inpnt  data  file  for  PREMOS.  When  the  PWL  analysis  is  needed,  the 
word  *pwl*  is  added  in  the  file  and  the  prograa  does  PWL  analysis  for 
the  DC  equilibrium  point. 


PLA  finite-state  aaohine  implementing  the  light  controller 

*subcircuit  model  card 

model  inv  nor2  (51  lOf  lOOf) 

model  nor3  andoi(5  5  1  lOf  lOf  lOf  lOOf  0  3) 

model  nor4  andoi(5  5  1  lOf  lOf  lOf  lOOf  0  4) 

model  notrl  trans(5  1  2  lOf  lOOf  lOf  50f  1  1) 

model  notr2  trans(5  1  2  lOf  lOOf  lOf  50f  2  1) 

model  notr4  trans(5  1  2  lOf  lOOf  lOf  50f  4  1) 

model  notr 5  trans(5  1  2  lOf  lOOf  lOf  50f  5  1) 

model  clkl  source  (4  1  lOn  5n  lOn  5n) 

model  clk2  source  (5  0  5n  5n  5n  5n) 

*  AND  plane 

xl  11  17  19  1  nor3 
x2  13  17  19  2  nor3 
x3  12  14  17  19  3  nor4 
x4  15  18  19  4  nor3 
x5  16  18  19  5  nor3 
x6  12  13  18  20  6  nor4 
x7  11  18  20  7  nor3 
x8  14  18  20  8  nor3 
x9  15  17  20  9  nor3 
xlO  16  17  20  10  nor3 

*  OR  plane 

xll  5  6  7  8  9  21  56  28  notr5 

xl2  3  4  5  6  22  56  29  notr4 

xl3  3  5  7  8  10  23  56  30  notr5 

xl4  6  7  8  9  10  24  56  31  notr5 

xl5  4  5  25  56  32  notr2 
xl6  1  2  3  4  5  26  56  33  notr5 
xl7  9  10  27  56  34  notr2 

*  output  registers 
xl8  28  35  55  49  notrl 
xl9  29  36  55  48  notrl 
x20  30  30  37  inv 

x21  31  31  38  inv 
x22  32  32  39  inv 
x23  33  33  40  inv 


x24  34  34  41  inv 
*  input  buffers 


x25 

57 

42 

55 

45  notrl 

x26 

58 

43 

55 

46  notrl 

x27 

59 

44 

55 

47  notrl 

•  input  registers 

x28 

45 

45 

50 

inv 

x29 

46 

46 

51 

inv 

x30 

47 

47 

52 

inv 

x31 

48 

48 

53 

inv 

x32 

49 

49 

54 

inv 

x33 

50 

50 

11 

inv 

x34 

45 

45 

12 

inv 

x35 

51 

51 

13 

inv 

x36 

46 

46 

14 

inv 

x37 

52 

52 

15 

inv 

x38 

47 

47 

16 

inv 

x40 

53 

53 

17 

inv 

x41 

48 

48 

18 

inv 

x42 

54 

54 

19 

inv 

x43 

49 

49 

20 

inv 

* input  sources 

▼si  55  0  clkl  01000100010001 
▼s2  56  0  clkl  00010001000100 
▼sO  57  0  clk2  1111000000111111 
vbO  58  0  clk2  1111000000000000 
▼cO  59  0  clk2  1111111111111111 


•snslysis  requests 
opt  1  1  3  1  1  1  1 
oontl  000 

preset  (35.5)  (36.5)  (55.1)  (56,1) 
tiue  In  In 
dc 
pwl 

plot  55  56  42  43  44  35  36 
plot  37  38  39  40  41  9  10 
plot  12345678 
•end  55  56  42  43  44  35  36 
send  7  9  37  38  39  40  41 
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*1  =  Vbe  9be+  'be 
*2  =Vbc*  9bc  +  'be 


Fig.  3(b)  PWL  Equivalent  Circuit  for  Bipolar  Transistor 
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F^g.  5  The  Circuit  Consists  of  4  Nand  Gates 
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Fig.  6(b)  Only  One  Resistor  Inserted  between. the  Gates  ■+. 
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11  s  Vqs  9gS  +  ^GS 

12  =  Vqo*9go+  ^gd 
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Fig.  10(c)  PWL  Equivlaent  Circuit  for  MOS  Transistor 


Circuit  Diagram  of  a  PLA 


